library(parallel)
library(tidyverse)
library(data.table)
library(feather)
library(sandwich)
library(lmtest)
library(ggplot2)
library(stargazer)
library(plm)
library(xtable)
library(lfe)

###Note: some i/o and cleaning is repeated from "3-Employees Main" at the top here 
everyone <- data.table(read_feather("Data/NY_Voters_and_Employees.feather"))
everyone[, minsalary:= min(BaseSalary.2011,BaseSalary.2012, BaseSalary.2013, BaseSalary.2014, na.rm=T), by=list(employeekey, PayBasis)]
everyone <- everyone[!(minsalary <7.25),]
everyone[, MaxSalary:= max(BaseSalary.2011, BaseSalary.2012, BaseSalary.2013, BaseSalary.2014, na.rm=T), by=list(employeekey, PayBasis)]

##add the undermin variable
everyone[!is.na(BaseSalary.2014) | !is.na(BaseSalary.2015) | !is.na(BaseSalary.2016), UnderMin:= FALSE] #everyone working in this period
everyone[(!is.na(BaseSalary.2014) | !is.na(BaseSalary.2015) | !is.na(BaseSalary.2016)) & (BaseSalary.2014 < 8.00 | BaseSalary.2015 < 9.00 | BaseSalary.2016 == 9.00| BaseSalary.2017 < 11.00) , UnderMin:= TRUE]


fit_mod_16 <- function(df){
  long <- data.frame(employee=rep(df$'employeekey',2),
                     underMin = rep(df$'UnderMin', 2),  voted = c(df$'Voted16', df$'Voted12'),
                     year = c(rep(16, nrow(df)), rep(12, nrow(df))))

model <- felm(voted~underMin*factor(year)|0|0|employee, data=long)

  return(model)
}


#function to fit 2014 v 2010 panels for a given dataset
fit_mod_14 <- function(df){
  long <- data.frame(employee=rep(df$'employeekey',2),
                     underMin = rep(df$'UnderMin', 2),  voted = c(df$'Voted14', df$'Voted10'),
                     year = c(rep(14, nrow(df)), rep(10, nrow(df))))
  
  model <- felm(voted~underMin*factor(year)|0|0|employee, data=long)
  return(model)
}


#function to fit 2017 v 2013 panels for a given dataset
fit_mod_17 <- function(df){
  long <- data.frame(employee=rep(df$'employeekey',2),
                     underMin = rep(df$'UnderMin', 2),  voted = c(df$'Voted17', df$'Voted13'),
                     year = c(rep(17, nrow(df)), rep(13, nrow(df))))

  model <- felm(voted~underMin*factor(year)|0|0|employee, data=long)
  return(model)
}



make_stats <- function(wagerate){
    placeboincreaseresults <- as.data.frame(matrix(ncol=8, nrow=1))
    colnames(placeboincreaseresults) <- c("cutpoint", "est", "pval", "p05", "CIlower", "CIupper", "n", 'se')

  	placeboincreaseresults$cutpoint[1] <- cutpoint <- wagerate
  	wageplacebodata[, UnderMin:= NULL]
  	wageplacebodata[!is.na(BaseSalary.2011) | !is.na(BaseSalary.2012) | !is.na(BaseSalary.2013) | !is.na(BaseSalary.2014), UnderMin:= FALSE] #everyone working in this period
  	wageplacebodata[(!is.na(BaseSalary.2011) | !is.na(BaseSalary.2012) | !is.na(BaseSalary.2013) | !is.na(BaseSalary.2014)) & (BaseSalary.2011 < cutpoint | BaseSalary.2012 < cutpoint | BaseSalary.2013 < cutpoint | BaseSalary.2014 < cutpoint ) , UnderMin:= TRUE]
  	everyone_under <- fit_mod_16(wageplacebodata[(wageplacebodata$PayBasis == "per Hour" & wageplacebodata$MaxSalary < (cutpoint + 7) & wageplacebodata$MaxSalary > (cutpoint -1.75)),]) #basically move the window up bit by bit.
  	
  	placeboincreaseresults$est <- summary(everyone_under)$coef[4]
  	placeboincreaseresults$pval <- summary(everyone_under)$coef[16]
    placeboincreaseresults$se <- sqrt(plm::vcovHC(everyone_under,typ='HC0', cluster='group')[3,3])
  	placeboincreaseresults$p05 <- ifelse(summary(everyone_under)$coef[16] < .05, 1, 0)
  	placeboincreaseresults$CIlower<- confint(everyone_under)[4]
  	placeboincreaseresults$CIupper <- confint(everyone_under)[8]
  	placeboincreaseresults$n <- summary(everyone_under)$df[2]
    placeboincreaseresults
} 


wageplacebodata <- everyone[!(everyone$UnderMin==T) & (!is.na(BaseSalary.2011) | !is.na(BaseSalary.2012) | !is.na(BaseSalary.2013) | !is.na(BaseSalary.2014)),]; dim(wageplacebodata); dim(everyone) #keep the same years as in the main UnderMin setup
placebocutpoints <- seq(from = 11, to= 40, length.out=1000)
placebocutpoints <- expand.grid(placebocutpoints, year=c(2016, 2014, 2017))


make_stats <- function(row){
    wagerate <- row[[1]]
    year <- row[[2]]
    placeboincreaseresults <- as.data.frame(matrix(ncol=9, nrow=1))
    colnames(placeboincreaseresults) <- c("cutpoint", "est", "pval", "p05", "CIlower", "CIupper", "n", 'se', 'year')
  	placeboincreaseresults$cutpoint[1] <- cutpoint <- wagerate
  	wageplacebodata[, UnderMin:= NULL]
  	wageplacebodata[!is.na(BaseSalary.2011) | !is.na(BaseSalary.2012) | !is.na(BaseSalary.2013) | !is.na(BaseSalary.2014), UnderMin:= FALSE] #everyone working in this period
  	wageplacebodata[(!is.na(BaseSalary.2011) | !is.na(BaseSalary.2012) | !is.na(BaseSalary.2013) | !is.na(BaseSalary.2014)) & (BaseSalary.2011 < cutpoint | BaseSalary.2012 < cutpoint | BaseSalary.2013 < cutpoint | BaseSalary.2014 < cutpoint ) , UnderMin:= TRUE]
    if(year ==2016){
  	   everyone_under <- fit_mod_16(wageplacebodata[(wageplacebodata$PayBasis == "per Hour" & wageplacebodata$MaxSalary < (cutpoint + 7) & wageplacebodata$MaxSalary > (cutpoint -1.75)),]) #basically move the window up bit by bit.
     }
    if (year ==2014){
       everyone_under <- fit_mod_14(wageplacebodata[(wageplacebodata$PayBasis == "per Hour" & wageplacebodata$MaxSalary < (cutpoint + 7) & wageplacebodata$MaxSalary > (cutpoint -1.75)),]) #basically move the window up bit by bit.
    }
    if (year ==2017){
       everyone_under <- fit_mod_17(wageplacebodata[(wageplacebodata$PayBasis == "per Hour" & wageplacebodata$MaxSalary < (cutpoint + 7) & wageplacebodata$MaxSalary > (cutpoint -1.75)),]) #basically move the window up bit by bit.
    }

  	placeboincreaseresults$est <- summary(everyone_under)$coef[4]
  	placeboincreaseresults$pval <- summary(everyone_under)$coef[16]
    placeboincreaseresults$se <- sqrt(vcov(everyone_under)[4,4])
  	placeboincreaseresults$p05 <- ifelse(summary(everyone_under)$coef[16] < .05, 1, 0)
  	placeboincreaseresults$CIlower<- confint(everyone_under)[4]
  	placeboincreaseresults$CIupper <- confint(everyone_under)[8]
  	placeboincreaseresults$n <- summary(everyone_under)$df[2]
    placeboincreaseresults$year <- year
    placeboincreaseresults
}


cl <- makeForkCluster(nnodes=4)
out <- parApply(cl, placebocutpoints, 1, make_stats)
out <- rbindlist(out)
fwrite(out, 'Intermediate Files/placebo_cutpoints.csv')
